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Abstract 

We examine several features of Bose-Einstein condensation (BEC) in an ex- 
ternal harmonic potential well. In the thermodynamic limit, there is a phase 
transition to a spatial Bose-Einstein condensed state for dimension D> 2. 
The thermodynamic limit requires maintaining constant average density by 
weakening the potential while increasing the particle number iV to infinity, 
while of course in real experiments the potential is fixed and N stays finite. 
For such finite ideal harmonic systems we show that a BEC still occurs, 
although without a true phase transition, below a certain "pseudo-critical" 
temperature, even for D=l. We study the momentum-space condensate frac- 
tion and find that it vanishes as 1/ \/N in any number of dimensions in the 
thermodynamic limit. In D< 2 the lack of a momentum condensation is in ac- 
cord with the Hohenberg theorem, but must be reconciled with the existence 
of a spatial BEC in D= 2. For finite systems we derive the iV-dependence 
of the spatial and momentum condensate fractions and the transition tem- 
peratures, features that may be experimentally testable. We show that the 
iV-dependence of the 2D ideal-gas transition temperature for a finite system 
cannot persist in the interacting case because it violates a theorem due to 
Chester, Penrose, and Onsager. 
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I. Introduction 

Bose-Einstein condensation (BEC) has been observed recently in several 
laboratories^! using magnetic traps to confine and cool alkali atoms. Related 
experiments on lithium have been reported by third group!. The number of 
atoms involved ranged from a few thousand to a few million, in potential wells 
that were to a good approximation anisotropic harmonic oscillator potentials. 
Such a spatial BEC gives rise to several theoretical questions. 

First, did these recent experiments observe a true phase transition to the 
BEC state? The answer is, obviously, no, because a true phase transition, 
with nonanalytic thermodynamic functions, requires taking the number of 
particles and the volume to infinity while keeping the density constant. Of 
course, no real system ever has such properties, but in most homogeneous 
macroscopic systems the thermodynamic limit is a good approximation to 
the experimental situation in which boundaries seem to play a relatively 
unimportant role. In the magnetic traps, not only is the number of particles 
quite small, compared to the usual case, but the "boundary," formed by 
a quadratic potential well, extends literally throughout the whole system. 
In order to take the thermodynamic limit in such a system it is necessary 
to weaken the potential so that, as the number of particles increases, the 
average density remains constant. This is well-defined mathematically, but 
is of course physically unrealizable. On the other hand, taking the box size to 
infinity in the homogeneous case is also unrealized experimentally. One can 
argue that the situation there is not qualitatively different from a gas in a 
harmonic trap. For the ideal gas all that matters is the density of states, and 
the thermodynamic limit simply carries that to a continuum in each case. 

The above discussion leads one to ask if the experiments then observed a 
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real Bose-Einstein condensation. The answer seems surely positive that they 
did see a macroscopic number of particles occupying the lowest harmonic 
oscillator state. Moreover the transition occurs quite abruptly as temper- 
ature is lowered. This result is in accordance with the findings of several 
authors§^0 that in finite homogeneous Bose systems there is an accumula- 
tion point or "pseudo-critical" temperature where the increase in the chem- 
ical potential slows and the number of particles in the ground state begins 
increasing rapidly. Many different definitions of this accumulation point have 
been offered for an ideal gas - all of which approach the true phase transition 
point in the thermodynamic limit. 

Thirdly, can a pseudo-transition occur in cases where there is no real 
phase transition in the thermodynamic limit? It can. One can showEl that 
there is a real transition in the harmonic potential for an ideal gas in any 
dimension greater than or equal to two. However, although there is no real 
transition in one dimension (ID), there is a pseudo-transition0'0 that occurs 
at temperature that would go to zero as 1/ In N in the thermodynamic limit. 

The results quoted in the last paragraph bring up yet another question. 
How can there be a transition in 2D in the harmonic potential when the 
well-known Hohenberg theorem@ says that there is no BEC transition in 
that number of dimensions? What the Hohenberg theorem actually says is 
that there can be no BEC into the k = state where k is wave number. 
Despite the restriction to /c-states, this theorem would seem to be relevant to 
the case of a BEC into a harmonic oscillator ground state, because there is 
yet another theorem (CPO theorem), due to Chesterffl, based on a lemma 
of Penrose and OnsageriLZl, that notes that there can be no BEC into any 
single-particle state unless one occurs into a k = state. So it seems that 
the 2D ideal gas transition ought not to be allowed! This situation has 
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arisen before and the answer founds The CPO theorem requires that the 
density be finite everywhere. In the thermodynamic limit, the density of an 
ideal gas becomes infinite at the origin in the harmonic oscillator problem, 
which negates the validity of the CPO theorem. So there can be and is a 
BEC into the harmonic oscillator ground state in 2D in the thermodynamic 
limit. 

But the Hohenberg theorem does not depend on the finiteness of the den- 
sity for its validity. So it must still be valid to say that, in the thermodynamic 
limit, there is no BEC into the k = state for the 2D oscillator problem, 
while there is one for the lowest oscillator state; that is, there is a spatial 
condensation but not a momentum condensation. We see then that we must 
define two condensate numbers: no, the number of particles in the lowest 
oscillator state, and fo, the number of particles in the zero-momentum state. 
If there are N particles in the system, then for the 2D ideal gas at non- 
zero temperature in the thermodynamic limit, fo/N approaches zero, while 
n /N remains a finite fraction. Basically these two quantities are Fourier 
transforms of one another, and, despite the presence of the harmonic poten- 
tial, it is still possible to have particles in the lowest k state. This peculiar 
relationship between them has not been noted previously. 

Our calculations are for the ideal gas. In the thermodynamic limit the 
density becomes infinite at the origin in a harmonic potential. If a hard-core 
interparticle repulsion were included, the density would not become infinite, 
the CPO theorem would apply and the 2D BEC in a harmonic potential 
would disappear in the thermodynamic limit. It is easy to extend the ideal 
gas calculations to include a mean-field interaction!!!!, but to the author's 
knowledge no one has yet done the 2D calculation for a system with hard- 
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core interactions in a harmonic potentials Because the actual experiments 
are done with just a few particles, the ideal gas should become a fair ap- 
proximation to the real experiments, and the density divergence of the ideal 
gas should never become an issue. Thus our calculations of no, fo, and 
the transition temperature might be relevant to experiment.-'' However, the 
CPO theorem tells us that the pseudo-transition temperature would vanish 
were we to take the thermodynamic limit. Thus the iV-dependence of the 
pseudo-transition temperature must be different from that of the ideal gas, 
meaning it must depend more weakly on particle number than for the ideal 
gas. This is in contrast to the 3D case where the ideal gas transition temper- 
ature seems to describe experiment quite wellH. A full hard-core calculation 
of the transition temperature in 2D would of course be useful. 

In Sec. II, we review the treatment of the Bose condensation in the ther- 
modynamic limit. We look at the finite case in more detail in Sec. III. Sec. IV 
examines the momentum distribution and its relation to the distribution in 
the harmonic states. Further discussion occurs in Sec. V. An Appendix gives 
some mathematical details and makes connection with previous work. 

II. Harmonic Bose Systems in the Thermodynamic Limit 

Bose-Einstein condensation in an external harmonic potential has been 
considered previously0~0'El~ii'ii~ii. We re-examine and extend the analysis 
here. Consider the two-dimensional system of N noninteracting Bosons with 
the particles contained in an isotropic two-dimensional harmonic potential 

TA 3D calculation of the transition temperature for a gas in a harmonic potential mth 
hard core interactions has been done recently by path-integral Monte Carlo methodstj. 

°^ As pointed out recently in Ref. [l3[ experiments on systems with reduced dimensionlity 
(ID or 2D) are entirely possible in magnetic traps. 
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U{r) = \uJ r -] 2 (1) 



given by 

" M = 2- V«7 

where r 2 = x 2 + y 2 and is a range parameter for the potential. It is simple 
to generalize to any number of dimensions greater than two. The energy 
levels are given 

E mx ,m y = hu(m x + m y + l) (2) 
with m x , m y = 0, 1, 2, . . . The angular frequency is 



(3) 



' R 2 m 

where m is the mass of a particle. In the grand ensemble at temperature T, 
the total number of particles N satisfies the equation 

OO -1 

iV= V f4) 

m x ,m y =0 

where /i is the chemical potential, (3 = 1/ksT, and ks is the Boltzmann 
constant. 

If one wishes to take the thermodynamic limit, then one must increase 
the "volume" of the system while increasing the number of particles with 
the average density kept constant. The average density is proportional to 
p = N/R 2 , where R is the range parameter in the harmonic potential of 
Eq. (H|). To take the thermodynamic limit requires that R increase, that 
is, that the potential weaken, as N increases. Such a situation has been 
considered previouslyB'lllS, but might seem peculiar when considering real 
experiments in a potential fixed by external magnets. We follow this pro- 
cedure, because it is the only way to reach the true nonanalytic behavior 
characteristic of phase transitions. One does not feel so uncomfortable with 
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the thermodynamic limit in the case of a homogeneous particles-in-a-box sit- 
uation because there it seems that all that the thermodynamic limit does is 
to remove boundary effects, which are mostly negligible anyway in the real 
bulk experiments. Making such a qualitative distinction between the two 
cases seems unwarranted, however. 

If p is kept constant we require, in 2D, that 

II 

uj oc — oc — = (5) 
In D dimensions, p oc and u oc N l /D . 



Introduce 



( 2 ) _ h u 



and 

a = —(3p + Uuj. (7) 

Then 



(2) 



k B T Ty/~N 

and Eq. (|J) can be written as 



(8) 



N= E "^5 " " (9) 

mx,rriv e ivw ( - mx+m y^ +a _ i 

To simplify, one can change the sums over m x and m y to sums over 
p = m x + m y and I = m x . Then, since I no longer appears in the summand, 
one can do that sum to give 

N = E J • (10) 

p=0 Q p+a i 

One can use this formula to do numerical calculations of a for given finite 
N. We will discuss this procedure in the next section. Here we want to take 
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the thermodynamic limit and it is useful to separate off the ground state by 
writing 

00 p+l 

n = ™ + E — 

p=1 e T7WP +a - 1 

^ P + 2 
= n o + l^— • (H) 

p=0 e T7WP+ a ' - 1 

where 

"o = ^y- (12) 

r (2) 

In the last form we have taken p — > p + 1 and a = a + T ^ to reset the sum 
from p = 0. 

When iV gets large, the states become very closely spaced and we can 
replace the sum by an integral (see Appendix) to good approximation to 
become 

N = n + dp^f (13) 

J T n _ , . / 



gTVJV — 1 



Changing variables to u — pT^ fy/NT, we find 

/ T \ 2 r°° u / T \ r°° 1 

N = n + N — TTyT / du — — h 2VN / du — — 

VT (2) / ^0 e«+ Q ' - 1 VT (2) / -/o e u+a ' - 1 

T \ 2 _ , „ / T 



where 



cr— 1 00 g — ap 



V ' r(or) Jo e u+a - 1 ^ 



(15) 

pu 

are the Bose integra isi; r( a) is the Gamma function. The term in F2 is of 
order iV while the term in F± is of order viV and can be neglected in the 
thermodynamic limit. We will however use it in the next section to estimate 
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the error made by the thermodynamic limit form when N is finite. Since F 2 
has a finite limit for small a, there is a BEC for temperatures lower than 
some critical temperature. F 2 (a) behaves for small a a s§ 

F 2 (a) =C(2)-a(l- In «) + ... (16) 

where ((a) = J2%Li 1/V 7 is the Riemann ^-function. We have C(2) = n 2 /6. 
The 2D Bose condensation temperatures is given by the equation 

N = N^F 2 (0). (17) 



or 



Tf> =T w ^«0.78lf ) . (18) 



The occupation of the lowest oscillator state is, for T < T^ 2 \ 

%- l -{wf 20 (19) 

and zero above T^ 2 \ Unlike the 2D homogeneous system the 2D oscillator 
has a BEC phase transition in the thermodynamic limit. 

Recently the authors of Ref. |29] have claimed that the harmonic potential 
is fundamentally different from a homogeneous system and that a true phase 
transition cannot occur in the harmonic potential. We see that this is not 
correct, although the physical conditions for its occurrence (weakening the 
potential) may seem a bit strange. However, the same arguments could be 
applied to the system of particles in a box. Real experiments do not occur 
in an infinitely large box and experimentalists do not observe actual phase 
transitions with truly discontinuous functions. Further, all the particles know 
is the energy levels and their spacing; in each case the limiting process is just 
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changing the spacing. We feel that the two cases are not fundamentally 
different. 

In 3D we have 

OO 1 

N= E (20) 

m x ,m y ,m z =0 — Q T7^( m x+m y +m z )+a 



(3) 

where T is given by 



T (3) = l\ % 1/3 ( 21 ) 



with p = A^/i? 3 and 

a = -fin + 3T (3) /(2TjV 1/3 ). (22) 

One can again reduce the sum of Eq. (^) to one variable over p = m x + 
m y + m z to give to 

N = no + t^k±L (23) 

with a' = a + To/fW 1 / 3 ). We have again separated off the no term and have 
reset the sum from p = 0. Changing the sum to an integral (see Appendix) 
gives 

N = n + N (^y)W) + (^)) 2 ^(«0 + 3^ 3 F^a') 

(24) 

In the thermodynamic limit the terms in F 2 and F% are negligible, as is the 
difference between a' and a, and the condensate fraction is given by 

N " It c (3) 



"" 1-fi) 3D (20) 



where = T (3) /(C(3)) 1/3 . 
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Eq. (|24| ) is very similar to one given recently in Ref. |2l, but differs in the 
coefficients of the various terms and in the argument being a 1 rather than a. 
We show in the Appendix that the two equations are the same. Our equation 



also has the same form as that appearing in Ref. |27. 

It is similarly easy to show that the D-dimensional harmonic oscillator 
system has a BEC, for D> 2, at temperature 



T (P) = T <P)£( D yi/D ( 26 ) 



with Tq D ^ given by 



T (D) = ^\-P 1/D (27) 
and p = N/R D . In D-dimensions (D> 2) the condensate fraction is given by 

— t = 1 — f — tttt | D dimensions. (28) 

N Vt c (d) J 

As oo, — > Tq D ' ) and uq/N becomes a step function in T. 

The ID system is a special case for which there is no condensation in 
the thermodynamic limit. The same procedure as used above leads, in the 
continuum limit, to 

N = N (^m) f ^7x^ = N (^m) *H«)- (29) 



/ Jo e u+a - 1 \T (1) 



with T (1) given by Eq. (|27|) with D= 1. -Fi(a) does not approach a finite 
limit as a — > 0, but is given exactly by 

Fx (a) = -ln(l -e- a ) (30) 



which approaches — In a as a — > 0. However, Refs. [12] and [L3] have pointed 
out that, in the finite system, the ID system does have a BEC below a certain 
pseudo-critical temperature. We will discuss this case below. 
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III. Finite Systems 

To consider finite harmonic systemsE'SB^i! most easily in, say, 2D, 
we can just compute a via Eq. ( |1T)D by iteration and then compute the num- 
ber of particles in the lowest state from Eq. fll2|) . Fig. 1 illustrates the 2D 
situation. The dotted and dashed lines are the exact results for N = 10 to 
10 4 with the thermodynamic limit (Eq. (|19|)) shown as the solid line. 

One notes that the exact result is smaller than the value given by the 
infinite N limit. A better approximation than Eq. (|I9| ) is the form of Eq. (|T4]) ; 
the F\ term in Eq. (|14]) aids in giving an estimate of the difference due to 
finite N. In 3D a better analytic approximation is given by Eq. (H)Bi-§. 
Since, for small temperature, a = 0(1/N), the F a behave as 

Fi(a) -> -ln(a) 

F 2 {a) -+ C(2) - a(l - ln(a)) (31) 
and we can show that the condensate fraction is roughly 

^~l-f^V-— ^ (32) 
N [tPJ T {2) VN 1 ] 

The correction is negative so that the value of no valid for finite N is smaller 
than the infinite limit value, at least when a is of order 1/N as we see in the 
figure. 

For a finite harmonic system it is not particularly appropriate to express 
the pseudo-transition temperature in terms of the density. The potential is 
fixed and scaling it in terms of iV with fixed density makes little sense. We 
now express it instead in terms of a fixed frequency. The effective or pseudo- 
transition temperature is only a bit less (less by terms of order In N/y/N) 
than the infinite-limit transition temperature and can be taken as roughly 
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the same. By Eqs. (|6]) and (|T^) the pseudo-transition temperature is given 
by 

A similar analysis for three dimensions leads to a pseudo-transition tem- 
perature a bit less than 

/jM" 3 ^ (34) 



A pseudo-transition in one-dimension is a special case!30 because there 
is no real phase transition in ID. Pseudo-transition temperatures can be 
defined in many ways. Ref. |7| defines six different ways in the homogeneous 
case that all lead to the same T c for infinite N in cases where there is a real 
phase transition. We compare just two approaches. First we follow^ the 
analysis used in the higher dimensional cases: Write 

OO 1 

N = E— (35) 

p=0 e -TwP +a - 1 



where Tq 1 ^ is given by 



T„"» , T-y—Ji (36) 
k B \ m R 



and a = —/3fi + / (2TN). Separate off n as above to get 



OO 1 

N = n + Y,— (37) 

p=0 



in which 



a' = a + ^ (38) 
and the sum has been reset from p = 0. The continuum limit becomes 

N = n + N-^F 1 (a'). (39) 



(i) 
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The assumption that no is of order N implies, from Eq. ([12]), that a is of 
order 1/N. We have, for small a', 



iV = n -^ y ln(a / ) (40) 



Taking a' of order 1/N and assuming N large give 



N = n + -^N\nN (41) 



T 



or 



I* = l (42) 

iV" 1 T «/lniV ( } 



This equation has just the form shown in Eq. (28) with an A^-dependent tran- 
sition temperature. (But, of course, in the case of a finite N and fixed har- 
monic potential all pseudo-transition temperatures are iV-dependent.) Here 
we have 

TW = I^ = — — (43) 
In AT lnjVfc B 1 ] 

in which the first form is appropriate to the case where the density N/R is 
being kept constant (weakening potential) and shows the transition temper- 
ature decreasing to zero as the number of particles increases so that there is 
no true phase transition. However, the second form shows that for constant 
potential as in the real experiments the temperature for which a substantial 
BEC occurs actually gets larger as iV increases. 

An alternative derivationl'El comes easily from noting that in a BEC one 
can write 

N = n (a) + N'(a) (44) 

where N' is the number of particles not in the ground state. At a real BEC 
transition temperature a is very small, but n is also still microscopic. Thus 
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at the transition temperature, we have approximately 

iV'(O) = N (45) 



In the case of the ID system this would imply, from Eqs. fl4~0|) and (|38| ) that 
is given by 

T (i) / T (D \ 
-±£-ln 4° =1 (46) 

or to a good approximation for large N, 

T (i) 

T (i) = io_ (47) 
In iV 1 ' 

in agreement with the first estimate. 

In Fig. 2, we plot the exact quantum result n /N from Eq. (^) as a 
function of T/Tq . Such a graph is appropriate to the case of taking the 
thermodynamic limit for which the density (and therefore T ) is held con- 
stant. One sees that then the transition temperature gets smaller as N 
increases, indicating that there is no real transition in this system. Fig. 3, 
on the other hand, plots Uq/N as a function of TlnN/T^. This method 
of graphing!! includes all the N dependence that occurs in and is more 
appropriate for the fixed potential case. One sees the asymptotic approach 



to an exact straight line for large N as given by Eq. (j42| ) 



IV. The Hohenberg and CPO Theorems 

The Hohenberg theorem0 states that there is no condensation into the 
k = state in the infinite system in two dimensions or less. The CPO 
theoremEl~ll3 says that there can be no condensation into any single particle- 
state unless there is one into the k = state. However, the CPO theorem 
relies on the assumption that the density is everywhere finite. This is not 
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the case in a harmonic potential, where we will show that the density be- 
comes infinite at the origin in the thermodynamic limit. We first review the 
Hohenberg theorem as presented by Chester^!. 

Chester has shown how to apply the Hohenberg theorem to a finite, inho- 
mogeneous system. To quote Chester!!!: "... [F]or any system whatsoever, 
we can always ask about the mean number of particles with momentum zero. 
This follows from the general principles of quantum mechanics and we are 
in no way limited, in asking this type of question, to systems that are ho- 
mogeneous." We consider the following picture. Assume that the harmonic 
potential extends out to radius R, but that there is a hard wall at r = R. It is 
possible, of course, to find the energy eigenvalues for such a mixed potential 
system, but we can use the eigenvalues of the harmonic potential if we just 
assume that k B T is small enough that particles rarely ever actually reach the 
hard wall. That is, we assume that R is sufficiently large that 

k B T < V(R) (48) 

The Bogoliubov inequality is 



h A A ti\ > k B T\{[C,A])\ 

r ' V- ([[c,#],ct]) 



2 



(49) 



in which Chester takes 



and 



A = a a\ (50) 



C = J2 a l a Q+k (51) 

q 

where a\ is the creation operator for the state with wave number k. To 
quote Chester once again0: "Since we are allowed to use any set of complete 
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functions for our second quantization formalism we choose a continuum of 
plane waves." Direct calculation of the commutators leads to 

, , , 1 / , , , ^ (/o ~ fkfk B Tm 
foh + ^(/o + A) > (52) 

where =< a k dk > is the number of particles in the momentum state k. 
This inequality is summed over k values from some minimum ko to maximum 
k c . The first term on the left can have the sum changed to all k values, which 
just increases the inequality and gives Nfo. The second term becomes NfoT, 
where T is of order unity. The third term is of order TV and is negligible if /o 
is assumed to be of order N. On the right side we can drop as negligible 
compared to / . Thus we have 

And so in 2D, we have (assuming that dropping the T term does not seriously 
affect the inequality) 

N > hmkBTA x K 
h 2 k 

where A = tcR 2 is the area of the system. We take k oc l/R and k c oc I /a, 



where a = l/y/p = \J R 2 /N. Then our result for the 2D case is 

ill < (55^ 

N ~ ma 2 k B T In N V ; 

Thus the Hohenberg theorem states that, in the thermodynamic limit, 
the condensate fraction in the k — momentum state must go to zero at 
least as fast as 1/lniV for any non-zero temperature. We will see that it 
actually goes to zero faster than this, namely as 1/ \^N. 

The CPO theorem goes as follows: Let <7i(r, r') be the single-particle 
reduced density matrix in the position representation. Penrose and Onsager's 
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lemma states that if the quantity 

W = vJ dr J rfr >i(r,r')|, (56) 

(where V is the volume) is not proportional to N, then there is no BEC into 
any eigenstate of the single-particle density matrix. However, the number of 
particles in the zero momentum state is given by 

fo = J dr J dr' < k = 0|r >< r|(7i|r' >< r'|k = > 

= ^ J dv J dv'a^v'). (57) 

Chester proves that a Bose system has a positive-semidefinite single-particle 
density matrix, so that, if fo/N approaches zero as N — > oo, then W/N must 
also go to zero, eliminating the possibility of BEC into any single-particle 
state. 

The weak link in the theorem's proof for our purposes is that the Penrose- 
Onsager lemma requires the average density at each point be finite. That is 
not the case in the harmonic potential as we will show. Thus the fact that 
fo/N — ■> does not imply that n /N — > 0. In the thermodynamic limit it is 
possible, for a 2D harmonic potential, to have a condensation into the lowest 
harmonic oscillator state while still satisfying the Hohenberg theorem that 
there be no BEC into the zero momentum state in the thermodynamic limit. 
We will indeed show that fo/N approaches zero as l/y/N. 

A direct computation of fk is straightforward. We have 

fk = ^2 < k| m >< m l°"il m >< m |k > (58) 

m 

where the \m > represent the eigenstates of Ui, in this case the harmonic 
oscillator eigenstates. The states < k|m > are the Fourier transforms of the 
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harmonic oscillator states. These latter functions are, just like the states in 
configuration space, Hermite functions times a Gaussian. In ID 



< k\m > = 



in which L is a box length, 



and 



moo 

IT' 



c„ 



V \/22 m m\ 

is the usual wave function normalization factor in configuration space. 
The density matrix is given in 2D by 

1 

v mx m y =< m x m y \a 1 \m x m y >- e/3[Mma;+m!/+1) _^ _ 1 
The result of the indicated calculations is 



(59) 

(60) 
(61) 



(62) 



A = ^e-(^) £ o mxmy ClCl\H m ^ x )H my ^ y )\ (63) 

? m x ,m y =0 

where Ki = and A is the area. Using the relations derived in Sec. II we 
can rewrite the last equation as 

4 h 2 



fk — 



Nk B T Q ma 

oo 

X £ °: 



mx , my =o " y 2^m x \m y \ 



H mx (K x )H my (K y ) (64) 



where a = y^l/p. 

The value of the zero momentum term is 
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We can find the density at position r by a similar analysis: 



oo 

e2„2 



p(r) =< r|ai|r >= e-« p £ a mxmy C 2 m C 2 my [H mx ^x)H my ^y)\ (66) 

mi,m,,=0 

Because of the relation between harmonic oscillator wave functions and 
their Fourier transforms, we see that there is a close relation between f and 
p(0). This is given by 

/° = 4 ^(^rV) p(0). (67) 

We can find WKB approximations for fk and p(r) by using the Boltzmann 
distribution function for the number of particles of wave number k at position 
r, given in 2D by 

fr B (r) = ~ A ^ l 7 - (68) 



If we integrate over r and sum over k we should get the number of particles in 
the system. However, in 2D we actually get just the F 2 (a) term in Eq. (fH]). 
This is N', the number of particles in excited states. We need to treat the 
condensate separately. 

If we integrate f^ KB ( r ) just over r we get f^ KB i the wave- number dis- 
tribution, and if we just sum over k we get p(r), the local density. The results 
of these integrations are 

^ = - J^T o)2 ln ( X " ^M-Pm 2 /2m - a}) (69) 

and 

p(r) = -^^ln(l-exp[-/3rj(r) -a]) (70) 

where a = —fi\x. 



21 



These formulas are valid only above any condensation or pseudo-critical 
temperature. We could add ^-function terms to account for the hq term, but 
we do not bother here. Were we to take k = and consider low temperatures 
such that a oc 1/N, the first of these functions would depend on — In a oc 
lniV. A similar behavior would occur for the second function for r = 0. 
We will see next, however, that the true quantum behavior gives a stronger 
divergence with N. 

In Fig. 4 we show plots of the quantum as a function of k as computed 
from Eq. (0). The lower temperature 0.5To is below the transition tempera- 
ture T c = 0.78To and the higher temperature is above T c . Also shown are the 
WKB approximations at the same temperatures from Eq. ( |69"D and Eq. (|70|) . 
We see that the WKB approximation is quite good if we are above T c , but 
very poor below as expected. 

One important feature of fk shown is that fo is continuous with for 
k 7^ 0. This is quite different from homogeneous Bose condensed systems 
for which the k = term is much larger that the k ^ terms. The pseudo- 
condensate in momentum space (terms not included in the WKB expression) 
is spread over many k-values. However, note that the /c-scale factor deter- 
mining the width of the spread in fc-space of the functions in the Fig. 4 
is 



which goes to zero as iV — > oo in the thermodynamic limit. 

In Fig. 5, we show the k = momentum occupation number as a function 
of temperature for iV = 1000. Clearly there is a pseudo-condensation at the 
same temperature T c = 0.78To as the spatial condensation into the oscillator 
ground state. 




(71) 
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The ^-dependence of fo is shown in Fig. 6. For the two lowest tem- 
peratures the behavior is found numerically to be fo oc y/N . The highest 
temperature shown is just below the transition temperature and has not yet 
settled into its asymptotic form. With this result we will have 

|oc^^0 (72) 

as iV — > oo. This supports our claim that the k — condensate fractions 
vanishes in the thermodynamic limit. 

We can verify this numerical result analytically. The first term in the 
sum in Eq. (|65|) is just <7oo = no so that the first term is of order y/N when 
no is of order N. The other terms in the series are positive so that the whole 
sum is at least of order y/N. The quantity in square brackets in Eq. (B^) is 
always less than 1. This implies that 

f ^ 4 ^ 1 V- f7 oy 

ma l£ B l Q y/N ~~ 



Since the sum on the right is just N, we have 



Ah 2 ,— 

fo < —nr^^ (74) 
ma^KBJ-o 



This result implies that fo/N approaches zero with iV as 1/y/N as our numer- 
ical treatment showed. Actually we can do somewhat better in our analysis. 
The quantity in square brackets in Eq. ( |65|) can be shown (via the Sterling's 
approximation) to obey 

m -- m v- 1 „ 2 (75) 

2 m -2 m -[(m x /2)\(m y /2)\] 2 \ n^m x m y 1 ; 

for large m x and m y . If we separate off the m x ,m y = term, use Eq. ([75|) 

for all m x ,m y , and then replace the resulting sum by an integral, we find 

fo^^^F^m + fr B (76) 
ma 2 k B T y/N 
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where f^ KB is gotten from Eq. ( |55D . The WKB term in Eq. (|76| ) diverges as 
In N. In this approximation the y/~N behavior comes from just the term in 
n /VN. We see from Eq. ([7§), that for low T, f oc y/N/a 2 T oc yfNp. Thus 
/o/iV goes to zero in the thermodynamic limit as 1/yN. However, when the 
potential remains fixed, so that the density increases with N, then /o oc N. 

Because of the relation Eq. (|67|), we see that the density at the origin also 
diverges as yN in the thermodynamic limit. 

It is possible to extend our discussion to arbitrary dimension. We find that 
in any dimension fo/N tends to zero as 1/yN so that in the thermodynamic 
limit there is no BEC into the k = state. While the Hohenberg theorem 
forbids a k = condensate only for for D< 2, there is, on the other hand, 
there is no theorem that requires a BEC into k = for D> 2. 

We should again point out that our results apply only to the ideal gas. 
For a real gas with hard-core interactions there can be no divergence in the 
density. If the particles have a hard-core interaction, the density will stay 
finite everywhere and the connection between the CPO and Hohenberg the- 
orem will be re-established. This implies that the 2D BEC into the lowest 
oscillator state will disappear for interacting particles in the thermodynamic 
limit. Thus the CPO theorem implies that for interacting finite 2D sys- 
tems the dependence of the transition temperature on particle number must 
change from oc yfp oc \fN as shown for an ideal gas to something more 
weakly dependent on N. An example would be the form T c oc y/N / In N 
analogous to the ID gas, but the determination of the exact form, of course, 
requires explicit hard-core calculations. Curiously, the 3D, and even the ID, 
dependences are not required to change. 
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V. Discussion 

We have examined several features of BEC in harmonic potentials for 
the ideal gas. We have shown that, in the thermodynamic limit, there is no 
condensation into the zero-momentum state in any number of dimensions. 
This result occurs in spite of the fact that there is a BEC into the lowest 
oscillator state for dimension D> 2. 

For 2D there is a peculiarity in finding a condensate in the lowest oscillator 
state in the thermodynamic limit given the usual connection, via the CPO 
theorem, between the k — condensation and the condensation into any 
single-particle state. The fact that the Hohenberg theorem requires that 
there be no real condensation via a phase transition into the k — state in 2D 
would normally rule out the possibility of BEC into any other state. However, 
we have seen that the density diverges at the origin in the oscillator case 
(in the thermodynamic limit), which invalidates the CPO proof. Thus while 
there is no BEC into the k = momentum state in 2D in the thermodymanic 
limit, there is one into the m = oscillator state. 

The thermodynamic limit is an approximation to reality that might seem 
less reasonable in a harmonic trap than in homogeneous systems. The num- 
bers of particles in the initial real experiments were quite small - only a 
few thousand to a few million. Further, the mode by which one must take 
the thermodynamic limit, weakening the potential while increasing the par- 
ticle number in such a way that the average density remains constant, might 
seem qualitatively different from enlarging the size of a box of particles at 
constant density. In the box case it would seem that the enlargement is 
merely a change in the boundary conditions, which affects the minority of 
particles at the walls. On the other hand, varying the size of the experimen- 
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tal harmonic trap affects every particle. However, looked in another way we 
can say that all that the particles in either case know is the energy levels and 
their spacing. In both the box and the harmonic trap the thermodynamic 
limit treats the transition to the case where the spacing becomes infinites- 
imally small. The process is thus completely analogous in the two cases, 
simply with differing densities of states. 

The experimental situation, whether it be with a box or a harmonic mag- 
netic trap, usually has a fixed potential (which of course can be altered) with 
possibly varying numbers of particles; obviously the thermodynamic limit is 
never approached. Thus it becomes more useful to look at finite systems and 
how their properties depend on the number of particles N. Although finite 
systems do not have actual phase transitions, they have pseudo-condensations 
in which a substantial fraction of the N particles fall into the lowest oscillator 
state below some temperature whose value depends on N. Real experiments 
involve such pseudo-transitions and we have shown that they can occur in 
any number of dimensions, even in ID. Experimentalists using magnetic traps 
may have the capability of manipulating the traps to simulate 2D or even 

idJI 

The finiteness of the systems in real experiments may be precisely what 
makes our results relevant. Because there are so few particles in the real sys- 
tems, they do indeed behave much like ideal gases. Our finite-size ideal gas 
model has a density that is well-behaved at the origin without the unrealistic 
divergent behavior shown in the thermodynamic limit. In our finite-sized 
model, we can investigate how such systems behave as iV varies. There are 
pseudo-condensations into both the lowest oscillator state and into the k = 
state, and the relevant variables have predicted dependences on N. The 
ID case is interesting because of the behavior of the transition temperature; 



26 



the case of 2D is particularly interesting from the point of view of the CPO 
theorem. It is possible that these N dependences could be measured ex- 
perimentally and the predictions demonstrated. What needs to be done in 
further research is to investigate the effect of hard-cores to check that, in the 
low density limit, the results we find persist. 

The CPO theorem can predict that the persistence of ideal gas behav- 
ior cannot be complete. If the particles have a hard-core interaction, the 
density will stay finite everywhere even in the thermodynamic limit and the 
connection between the CPO and Hohenberg theorem will be re-established. 
This implies that the 2D BEC into the lowest oscillator state will disappear 
for interacting particles in the thermodynamic limit. The theorem thus im- 
plies that for interacting finite 2D systems the dependence of on particle 
number must change from that shown by an ideal gas in such a way that it 
can vanish in the thermodynamic limit. 

We note that is is possible that the prediction for f Q could be checked 
rather directly in experiments using a common method of measuring the 
condensate. The magnetic trap is removed and the particle cloud expands. 
Its "shadow" is then examined after it has blown up to a size large enough 
to be easily visible. The spatial distribution after removal of the trap is 
determined by the velocity distribution. The particles staying at the center 
of the cloud are those with zero momentum. 

Appendix 



Eqs. (|14|) and ( 24]) give corrections to the large-iV limit of the relation 
giving no(T,N). It is useful to consider a more rigorous derivation than we 
have given above and to make connection with previous work on the subject. 
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In particular, Eq. (B3j) looks different from the equivalent formula of Ref. |23 
[Eq. (4) in the first paper of the reference], in having for example a I in the F 2 
term instead of a |. We will see however that our formula agrees completely 



with the results of Ref. |24. 

We use the Euler-Maclaurin summation formulail to change our sums 
to integrals. There is more than one form of this formula. Eq. (3.6.28) of 



Ref. 32 



gives 

n-l 



£/(*) = / f(k)dk--[f(0) + f(n)} + -[f(n)-f'(0)} 

- ^[f»-r(0)] + --- (77) 



Note that the sum starts at k — 1. Eq. (|20|) can be written as 



00 1^2 1 3. 



p=l C 1 



where 



(3) 







TiV 1 / 3 ' 



(79) 



Use of Eq. (|77D leads to 



N = n + ^F 3 (a) + ~F 2 (a) + ■■■ (80) 

exactly as in Ref. However, the corrections to the Euler-Maclaurin in- 
tegral form a divergent series. It is easy to see that the various derivative 
terms are in powers of Uq = l/[exp(a) — 1], which for low T is of order N. 

To find a correction series that converges we use a slightly different form 
of the Euler-Maclaurin formulai; Eq. (23.1.30) of Ref. [32] is 

£/(*) = / o f(k)dk + -[/(0) + /(n)] + -[f'(n) - f 



k=0 



1 [/» - no)] + • • • (si) 



720 
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This form, in which the sum starts at k = 0, can be used directly on Eq. (|23|). 
Using it we get 



N 



3 

+ t 



1 



2 e Q ' - 1 



1 

12 



5/2 



+ ■ 



(82) 



e a ' - 1 (e a ' - 1) 2 _ 
Because a' = a + 6 appears in the arguments instead of a the series of 
corrections is well-behaved as we now see. 

First, we consider the case when T > T c so that a is of order unity while 
b remains of order iV -1 / 3 . We can simply Taylor expand around a' = a. We 
make use of the easily proved relation 

d 



da 



FJa) = -F a Ja) 



(83) 



for a > 1. After a bit of algebra we find to the same order 
N 



n + -F 3 {a) + ~F 2 (a) + -F^a) 



31 1 2 

H nn H n n 

24 12 



s 31 A 



2 



(84) 



Since we are above T c , rio is of order unity; b is of order iV x / 3 so that the 
entire series in various powers of Uq and b is negligible. We have then 

1 



N 



1 „ , n 3 1 , , 
/3(«)+2/ 2 («) 



(85) 



in agreement with Refs. |24|. Ref. has given a similar analysis. 

Finally, below T c , a is of order 1/N, which is much smaller than b, so 
that a' ~ b and 

1 l/6 = 0(iV 1/3 ). (86) 



e a ' - 1 

The full Taylor expansion then gives 



3 1 



N = n + N-F 3 (0) + --F 2 (0) + O - 



fe 3 



2 6 2 ' 



(87) 
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The first three terms correspond to the formula of Ref. |24], which is therefore 
valid below the transition as well as above. 
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Figure Captions 



Figure 1: Spatial condensate fraction n /N versus temperature (T/T ) in 2D. 
The dotted and dashed lines are for finite N values, generated from Eq. fllPP; 
the solid curve is the thermodynamic limit, Eq. (|l9|) . To is given by Eq. @. 



Figure 2: Spatial condensate fraction uq/N versus T/Tq in ID from Eq. (|35|). 
This curve is appropriate for considering the thermodynamic limit for which 
the average density and thus To are held constant by weakening the potential 
as N increases. The pseudo-critical temperature then decreases with N as 



1/lniV. In ID T is given by Eq. (37) with D=l. 



Figure 3: Spatial condensate fraction uq/N versus (T/To)lniV in ID from 
Eq. (|35|). This plot removes the N dependence from the approximate linear 
expression, Eq. ([12|). As N increases the exact expressions approach the lin- 
ear limit. When the harmonic potential is fixed while iV increases, as in real 
experiments, Tq/IylN oc iV/lniV; that is, the pseudo-critical temperature 
increases with N. 



Figure 4: Reduced momentum distribution function versus reduced mo- 
mentum component k x /£ in 2D. Here f% = fk[ma 2 kBTo/h 2 } and we have 
taken k y = k x . The upper solid curve is the exact quantum results from 
Eq. (|64j) for a temperature above T c and the lower one for a temperature 
below T c . The dotted curves are the WKB approximations at the same tem- 
peratures as the solid curves. The WKB expression is a good approximation 
only if T > T c . 
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Figure 5: Reduced zero-momentum distribution /q versus T/Tq in 2D from 
Eq. @5|) in 2D for JV = 1000. Here f* k = f k [ma 2 k B T /h 2 ]. A pseudo- 
condensate in momentum space clearly occurs for T <T c k, 0.78T . 

Figure 6: Reduced zero-momentum distribution /q versus A" in 2D [Eq. (|65[) 1 
for several temperatures below T c . These are log-log plots, with slopes of the 
upper two curves numerically very close to 0.5, indicating that fo oc y/N. 
The lowest curve, for T just below T c , has not yet reached its asymptotic A" 
dependence. 
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